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SOURCE CODE 

/* procedure: dsm_rindex . sas revision: 01 */ 

/* authored_by: doug mcnair modified_by: _ */ 

/* created_date: 04-JAN-2003 modif ied_date : */ 

/* purpose: data transform for adjusting health service */ 

/* for variations in rurality distance or access */ 

/* environment: H:\cerdsm\NHS\test */ 
/* =========================================================== ^ ====lt/ 

/* */ 
/* Read in data file. Assign the name "pop" to catchment pop.. : */ 
/* Assign the name n d_nhs" to distance to care category to.be */ 
/* distance- index adjusted. */ 
/* */ 
/* Only the infile statement and the the input statement need. to */ 
/* be modified. . .; ..*/ 

/* */ 
/* This is a locally normal index and assigns a quantitative . */ 
/* measure of rurality to each participant in a study. Data req'd */ 
/* to compute the index are minimal; only two variables need be. */ 
/* measured. Analysis effort required to compute the index is *./ 
/* negligible, and the tradeoff between data collection, and data */ 
/* analysis is defensible. Data collection is generally more . */ 
/* expensive than data analysis so minimize that. The two required*/ 
/* variables, county population and distance to care location, are*/ 
/* optimally transformed to achieve normality and weighted for */. 
/* validity. Measure of departure from normality is automatically */ 
/* obtained while constructing the index. Positive score reflects */ 
/* rural loc or urban loc with compromised access. Negative score */ 
/* reflects an urban residence with good access relative to the */ 
/* group under study. A score of zero reflects average distance */ 
/* or accessibility of health services, rural or otherwise. */ 
/* Reliability and construct validity were examined using two: */ 
/* data sets. */ 
/* */ 
/* The choice of a distance measure to gauge the distance between */ 
/* two distributions, in this case Gaussians, is critical.. .*/ 
/* The Kullback-Leibler, Bhattacharyya and arc-cosine distances */ 
/* are all related to maximum likelihood, and reduce to the ■-.-*/ 
/* Mahalanobis distance when used to measure the distance between */ 
/* two Gaussian distributions with equal covariance matrices. */ 
/* The advantage of the arc cosine distance is that it has a */ 
/* closed form expression, implemented in IML. */ 

/* */ 
/*======= =======================m 

data distance; 

infile 'c:\0_cerdsrn\health_econ\rurality\cancerpop.dat' ; 
input case__id d_nhs pop; 
retain d_nhs pop; 
use=l; 

if d_nhs = 0 then d_nhs= . 5 ; 
if d__nhs =. then use=0; 
if pop = . then use=0; 
keep d_nhs pop use; 
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if use = 1 then output; 
proc means; 
proc freq; 

tables d_nhs; 

proc iml; 

use distance; reset noname; read all var _num_ into YY; 
step=.05; mino=step; maxo=4 . ; guess={l, l}; 
k=2; c={.5, -l}; n=nrow (YY) ; 

do i=l to k; YY[,i]=YY[,i]/max(YY[,i] ) ; end; 
j=(l:n)\- p=(j-.5)/n; pi=arcos ( -1) ; 
X=J(n,k,0) ; Xl=J(n,k,0) ; X2=J(n,k,0) ; 
eps=l/2**26; 



start deriv; 

tau=lam-l; omega= (ssq(tau) ) ; romega=sqrt (omega) ; IMP=I (k) ; 
if romega >= eps then do; 
tau=tau/romega; IMP=IMP-tau*tau * ; end; 
do j=l to k; 
y=YY [ , j ] ; 

if lam[j] =0 then do; 
t=log(y) ; 
tl=t#t/2; 

t2=2*tl#t/3; end; 
else do; 

t=sign (lam [ j ] ) *y##lam [ j ] ; 

tl=t#log(y); t2=tl#log (y) ; end; 
t=t-J(n / l / sum(t)/n) ; vO=ssq(t) ; xj =t/sqrt (vO) ; X[,j]=xj; 
tl=tl-J (11,1, sum (tl)/n) ; t2=t2-J(n,l,sum(t2)/n) ; 
v01=xj**tl; v02=xj**t2; 

xlj=(tl-xj*v01)/sqrt(v0) ; XI [ , j ] =xlj ; vll=xl j * *tl ; 
X2j=(t2-xj*v02-2*xlj*v01-xj*vll)/sqrt (vO) ; X2[,j]=x2j; 
end; 

z=X*c; s2=z**z/ (n-1) ; s=sqrt(s2); z=z/s; 
ii=rank(z); B=z; z[ii]=B; 

B=X; X[ii / ]=B; B=X1 ; Xl[ii,]=B; B=X2 ; X2[ii,]=B; 
G=sqrt (omega/s2) *diag(c) * (XI*- (Xl**z) *z7(n-l) ) ; 
Phi=probnorm(z) ; 

Phi = Phi#(eps <= Phi) # (Phi <= 1-eps) + J (n, 1, eps) # (Phi < eps) 

+ J (n,!, 1-eps) # (1-eps < Phi) ; 
logphi=log(Phi) ; logphil=log (1-Phi) ; 

ph=exp(-z#z/2)/sqrt (2*pi) ; logph=- (z#z+J (n, 1 , log (2*pi) ) ) /2; 
A2n=2 *p * *logphil-2*p** logphi -2 *sum( logphil) -n; 
al=-2*exp ( logph- logphi -logphil ) # (p- Phi ) ; 

a2=2*exp (2*logph-2 *logphi -2* logphil) # ( (p-Phi) # (p-Phi) +p# (1-p) ) ; 
a2=a2+2*exp (logph- logphi -logphil) #z# (p-Phi) ; 
dl=G*al; f=IMP*dl; 

d2=(omega/s) *diag(c) * (X2 x *al- (X2**z) * (z^al) /(n-1) ) ; 
d3=sqrt (omega/s2) *IMP*diag (c) *X2 * *z ; 

H=- (f *tau*+tau*f *) ; Hl=diag(d2) - (tau**dl) *I (k) - (dl*d3 * +d3*dl *) / ( 

H1=H1-G*G* * (z * *al) / (n-1) + (G# (J (k, 1, 1) *a2 * ) ) *G X ; 

H=H+IMP*H1*IMP; 

finish; 

start optimize; 
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delta=l; A2n0=n; 

do while. (delta > l.e-5); 

run deriv; 

call gsorth(T,u, lindep, I (2) -tau*tau* ) ; T=T[,1]; 

b=-inv(T**H*T) * (T**f ) ; 

tau=tau+T*b; 

tau=tau/sqrt (ssq(tau) ) ; 

1 am= 1 + 1 au * r omega ; 

delta=n*abs (A2n0-A2n) ; 

A2nO=A2n; 

end; 

finish; 



start search; 

filename out ' power. out'; 

file out; 

lam={l,l}; run deriv; 

put @2 romega 7.4+3 A2n 9.4+3 (lam[l] ) 7.4+3 (lam [2] ) 7.4+3; 
lam=guess- . 001; 

do romega=mino to maxo by step; 
omega=romega**2 ; 

tau=lam-l; tau=tau/sqrt (ssq(tau) ) ; 1 am=l+ romega *tau; 
run optimize; 

put @2 romega 7.4+3 A2n 9.4+3 (lam[l] ) 7.4+3 (lam [2] ) 7.4+3; 
end; 

closefile out; 
finish; 

lam={0.2, - .1}; run deriv; print lam A2n; 
lam={0 / - . 15} ; run deriv; print lam A2n; 
lam={0,-.2}; run deriv; print lam A2n; 
lam={l.,l.}; run deriv; print lam A2n; 
lam={0.,0.}; run deriv; print lam A2n; 
run search; 
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/* ==================== : ====== = = = ^ 

/* procedure: dsm_riskadj . sas revision: 01 */ 

/* authored_by: doug mcnair modif ied_by: */ 

/* created_date : 04-JAN-2003 modif ied_date : */ 

/* purpose: calculate risk-adjusted area rates for indicators*/ 

/* variations in rurality distance or access */ 

/* adj for age and dcat */ 

/* environment: H:\cerdsm\NHS\test */ 

/•""" —-———-*/ 

filename control 'h:\cerdsm\nhs\control.sas'; 



% include control; . 

%macro area; 

%if fcmsalevl eq 1 %then msa; 
%else hospstco; 
%mend ; 

%macro cty2msa; 

%if fcmsalevl eq 1 %then %do; 

label msa ='msa %area assigned from hospstco'; 
msa = input ( put (hospstco, area .) , 5.); 
%end; 
%mend ; 

%macro sum2msa; 

%if &msalevl eq 1 %then %do; 
proc summary nway; 

class msa popcat dcat; 
var pop; 

output out=qipop (drop=_type freq_) sum=; 

run; 
%end; 
%mend ; 

/* * 

/* add pop denominator * 
/* * 

proc sort data=in3 . &inf ilea3 . (keep=%area) out=%area nodupkey; 

by larea; 

run; 



data qipop; 

length hospstco 8 popcat 3 dcat $1 ; 

infile popfile3 missover; 

input hospstco popcat dcat pop; 

%cty2msa 

run; 

%sum2msa 
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proc sort data=qipop; 
by %area popcat dcat; 
run; 

data qipop; 
merge %area ( in=x) qipop ; 
by %area; 
if x; 
run; 

/* it/ 

/* select obs for each indicator from previous step in turn. */ 
/* each iteration passes vars that control covars for that indie: */ 
/* n - obs number from covar files ■ . */ 
/* also identifies format used to index covars . */ 
/* pq - indicator name without prefix *./ 
/* a - number of demographic categories (age and d_nhs) . */ 
/* s - determines whether age*d interactions are needed */ 
/* o - determines omitted age cat in format *•/ 
/* it/ 

%macro mod3 (n,pq, a, s, o) ; 

/* .„ it/ 

/* SET creates TEMPI, containing the dep var TPQ and indep vars . *./ 
/* used in the regression. append to enctr data one obs per %area *■/. 
/* and demographic grp * / 

/* ^ 

data tempi (keep=key %area t&pq. d agel-age%eval (2*&a . ) aidx pop 
merge var) ; 
length agel-age%eval (2*&a . ) aidx 8; 

set in3 . &inf ilea3 . (keep=key %area t&pq. dcat agecat age) 
qipop; 

/* 

/* retain only enctr/area/demographic grp needed for this indicator*/ 
/* 

n = &n. ; 

if %area eq . then delete; 

if pop > 0 and n not in (2,9) then t&pq. = 0; 

if agecat in ('1') or popcat in (1,2,3,4) then do; 
if n in (1,3,5,7,8,13,14,15,16) 
then t&pq. = . ; 

end; 

else if agecat in ('2') or popcat in (5,6,7,8) then do; 
if n in (4, 6, 9) 
then t&pq. = . ; 

end; 

else if agecat in ('3', '4') or 

popcat in (9,10,11,12,13,14,15,16,17,18) then do; 
if n in (4,6,9) 
then t&pq. = . ; 
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end; 

if t&pq. ne . ; 

z* _____ v 

/* initialize covar dummies to zero */ 
/* ir/ 



array agex{ %eval (2*&a . ) } agel-age%eval (2*&a . ) ; 
d = 0; 

do i - 1 to %eval(2*&a.) ; 

agex(i)=0; 
end; 



/* it/ 

/* set index for age array (AIDX) . categories vary by indicator and*/ 

/* are determined by format corresponding to N */ 

/* - . ± j 



if age > 0 then aidx = put (age, a&n . f mt . ) ; 
else if pop > 0 then do; 
aidx = popcat; 

if aidx = &o. then aidx = 0; 
end ; 

else aidx = 0 ; 

/* _ ie/ 

/* set appropriate covariate dummy to 1 for this enctr */ 
/* _ ^ 7 

if dcat in (&s.) then d = 1; 

if aidx > 0 then do; agex(aidx) = 1; agex (&a . +aidx) = 1 * d; end; 
mergevar = 1 ; 



run; 

proc sort data=templ ; 
by %area d aidx pop; 
run; 



/* r it/ 

/* subtract number of encounters from %area/demographic group pop */ 
/* and set weight of 1 for each obs from enctr data, and a */ 
/* weight of (pop -encounters) for the appended obs for each group. */ 
/* it/ 

data tempi (keep=key larea t&pq. d agel-age%eval (2*&a . ) wht mergevar); 

set tempi; 

by larea d aidx; 

retain went 0; 

if first. aidx then went = 0; 

if pop > 0 then wht = pop - went; 

else wht = 1; 

went + 1; 
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run; 

/* ie/ 

/* read overall mean file for indicator; F variables are filler- */ 
/* ic/ 

data templm(drop=n f 1-f %eval (40- (2*&a . ) ) ) ; 

length n f 1-f %eval (40- (2*&a. ) ) 3 md mage 1 - mage %eval (2*&a . ) 8; 
infile meansa missover lrecl=1000; 

input n md magel-mage%eval (2*&a . ) f 1-f %eval (40- (2*&a . ) ) ; 
mergevar = 1; 
if n = &n,; 

run; 

/* + i 

/* subtract overall mean of vars from the vars on each obs ■ -.*/' 
/* so that adj area rate is the rate for the average person */ 
/* _ ie/ 

data tempi (keep=key %area t&pq. d agel-age%eval (2*&a . ) wht) ; . 

merge tempi templm; 

by mergevar; 

array covar{%eval (1+ (2*&a. ) ) } d agel-age%eval (2*&a . ) ; 

array means{ %eval (1+ (2*&a . ) ) } md mage 1 - mage %eval (2*&a. ) ; 

do i = 1 to %eval (1+ (2*&a. ) ) ; 

covar(i) = covar(i) - means(i); 
end ; 

run; 

/* ie/ 

/* read regression coeffs for each covariate */ 
/* „-„*/. 

data temp2 (keep=intercep d agel-age%eval (2*&a . ) _type ^model_ 

_depvar_ t&pq.); 

length intercep d agel-age%eval (2*&a . ) f 1-f %eval (40- (2*&a . ) ) t&pq.. 8 

_type model depvar_ $8; 

infile covara missover lrecl=1000; 

input n intercep d agel-age%eval (2*&a . ) f 1-f %eval (40- (2*&a. ) ) 
_depvar model type rmse_ t&pq.; 

if n = &n. ; 

run; 

/* 4r/ 

/* calculate predicted indicator for each obs using actual covars */ 

/* and the estimated coeffs, on de-meaned sample, where %area : */ 

/* effects are removed by subtracting the predicted val from the */, 

/* actual on each obs, and calculating the mean difference for the */ 

/* %area. See Green WH. Econometric Analysis 2e. 1990 pp. 466-9. */ 

/* it/ 
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proc score data=templ score=temp2 type=parms out=temply ; 

var d agel-age%eval (2*&a. ) ; 
run ; 

data temply; 

set temply; 

ehat = t&pq. - mhat; 
one = 1 ; 

run ; 

proc summary data=temply ; 
class %area; 
var ehat one; 

output out=r&pq. mean (ehat) =r&pq. sum (one) =p&pq. ; 

weight wht; 

run; 

%mend ; 

%M0D3 (1,APQ01, 18, 'F' ,6) ; 
%M0D3 (2 , APQ02 , 4 , ' F' ,0) ; 
%M0D3 (3 , APQ03 , 18 , ' F' , 10) ; 
%M0D3 (4 # APQ04, 18, 'F' , 1) ; 
%M0D3 (5,APQ05, 18, 'F' , 12) ; 
%M0D3 ( 6 , APQO 6,18, 'F',1); 
%M0D3 (7,APQ07, 18, 'F' , 10) ; 
%M0D3 (S^PQOS, 18, 'F' ,13) ; 
%M0D3 ( 9 , APQO 9 , 2 , ' F ' , 0 ) ; 
%MOD3 (10 , APQ10 , 18 , ' F' , 1) ; 
%MOD3 (11 , APQ11 , 18 , ' F' ,5) ; 
%MOD3 (12,APQ12,18, 'F' ,6) ; 
%MOD3 (13 , APQ13 , 18 , ' F' ,11) ; 
%MOD3 (14,APQ14,18, 'F' ,8) ; 
%MOD3 (15 , APQ15 , 18 , ' F' , 7) ; 
%MOD3 (16 , APQ16 , 18 , 7 F' , 11) ; 



data riskadj ; 

length %area RAPQ01-RAPQ16 PAPQO 1 - PAPQ16 8; 

merge RAPQ01 (KEEP=%AREA RAPQ01 PAPQ01) 

RAPQ02 (KEEP=%AREA RAPQ02 PAPQO 2) 

RAPQ03 (KEEP=%AREA RAPQ03 PAPQO 3 ) 

RAPQ04 (KEEP=%AREA RAPQ04 PAPQ04) 

RAPQ05 (KEEP=%AREA RAPQ05 PAPQO 5) 

RAPQ06 (KEEP=%AREA RAPQ06 PAPQO 6) 

RAPQ07 (KEEP=%AREA RAPQ07 PAPQO 7) 

RAPQ08 (KEEP=%AREA RAPQ08 PAPQO 8) 

RAPQ09 (KEEP=%AREA RAPQ09 PAPQO 9) 

RAPQ10 (KEEP=%AREA RAPQ10 PAPQ10) 

RAPQ11 (KEEP=%AREA RAPQ11 PAPQ11) 

RAPQ12 (KEEP=%AREA RAPQ12 PAPQ12) 



/* 



merge %area adjusted rates for indicators 



*/ 
*/ 
*/ 
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RAPQ13 (KEEP=%AREA RAPQ13 PAPQ13) 

RAPQ14 (KEEP=%AREA RAPQ14 PAPQ14) 

RAPQ15 (KEEP=%AREA RAPQ15 PAPQ15) 

RAPQ16 (KEEP=%AREA RAPQ16 PAPQ16) 
%area; 



label 

RAPQOl = 'DIABETES SHORT TRM COMPLICATN (Risk Ad j ) ' 

RAPQ02 = ' PERFORATED APPENDIX (Risk Ad j ) ' 

RAPQ03 = 'DIABETES LONG TERM COMPLICATN (Risk Ad j ) ' 

RAPQ04 = ' PEDIATRIC ASTHMA (Risk Adj ) ' 

RAPQ05 = 'COPD (Risk Adj)' 

RAPQ06 = 'PEDIATRIC GASTROENTERITIS (Risk Adj)' 

RAPQ07 = 'HYPERTENSION (Risk Adj)' 

RAPQ08 = 'CONGESTIVE HEART FAILURE (Risk Adj)' 

RAPQ09 = ' LOW BIRTH WEIGHT (Risk Adj)' 

RAPQ10 = 'DEHYDRATION (Risk Adj)' 

RAPQ11 = 'BACTERIAL PNEUMONIA (Risk Adj)' 

RAPQ12 = 'URINARY INFECTION (Risk Adj)' 

RAPQ13 = 'NEW CANCER MGT DELAY>14D S/P (Risk Adj)' 

RAPQ14 = 'DIABETES UNCONTROLLED (Risk Adj)' 

RAPQ15 = 'ADULT ASTHMA (Risk Adj)' 

RAPQ16 = 'LOWER EXTREMITY AMPUTATION (Risk Adj)' 



label 

PAPQ01 = 'DIABETES SHORT TRM CMPLICTN 

PAPQ02 = 'PERFORATED APPENDIX 

PAPQ03 = 'DIABETES LONG TRM CMPLICTN 

PAPQ04 = 'PEDIATRIC ASTHMA 

PAPQ05 = 'COPD 

PAPQ06 = 'PEDIATRIC GASTROENTERITIS 

PAPQ07 = 'HYPERTENSION 

PAPQ08 = 'CONGESTIVE HEART FAILURE 

PAPQ09 = 'LOW BIRTH WEIGHT 

PAPQ10 = 'DEHYDRATION 

PAPQ11 = 'BACTERIAL PNEUMONIA 

PAPQ12 = 'URINARY INFECTION 

PAPQ13 = 'NEW CANCER MGT DELAY>14D S/P 

PAPQ14 a 'DIABETES UNCONTROLLED 

PAPQ15 = 'ADULT ASTHMA 

PAPQ16 = 'LOWER EXTREMITY AMPUTATION 



(Denominator) 
(Denominator) 
(Denominator) 
(Denominator) 
(Denominator) 
(Denominator) 
(Denominator) 
(Denominator) 
(Denominator) 
(Denominator) * 
(Denominator) 1 
(Denominator) ' 
(Denominator) 1 
(Denominator) 1 
(Denominator) ' 
(Denominator) ' 



array arryl{l6} rapq01--rapql6 ; 
do' i s i to 16; 

if arryl(i) < 0 then arryl(i) = 0; 

else if arryl(i) > 1 then arryl(i) 
end; 



= l; 



drop i; 



run; 



/* 

/* smoothed demeaned rates = area rate * (signal / signal + noise) 
/* 
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data smoothed (DROP=RAPQ01-RAPQ16) ; 
set riskadj ; 



% include msx; 

array arrya{l6} RAPQ01-RAPQ16 

array arryn{l6} PAPQ01-PAPQ16 

array arryx{l6} SAPQO 1 - SAPQ1 6 



do i = 1 to 16; 

arryx(i) = arrya(i) - arrya3(i); 
if arryn(i) ne 0 and arrya2(i) ne 0 
arryx(i) = arryx(i) * (arrya2(i) 
arryn(i) ) ) ) ; 

else arryx(i)=.; 

arryx(i) = arryx(i) + arrya3 (i) ; 
end; 



then 

/ (arrya2(i) 



+ (arryal(i) / 



label 

SAPQ01 = ' DIABETES SHORT TRM COMPLICATN (Smoothed)' 

SAPQO 2 = 'PERFORATED APPENDIX (Smoothed)' 

SAPQO 3 = 'DIABETES LONG TERM COMPLICATN (Smoothed)' 

SAPQO 4 = 'PEDIATRIC ASTHMA (Smoothed)' 

SAPQO 5 = 'COPD (Smoothed)' 

SAPQO 6 = 'PEDIATRIC GASTROENTERITIS (Smoothed)'' 

SAPQO 7 = 'HYPERTENSION (Smoothed)' 

SAPQO 8 = 'CONGESTIVE HEART FAILURE (Smoothed)' 

SAPQO 9 = ' LOW BIRTH WEIGHT (Smoothed)' 

SAPQ10 = 'DEHYDRATION (Smoothed)' 

SAPQ11 = 'BACTERIAL PNEUMONIA (Smoothed)' 

SAPQ12 = 'URINARY INFECTION (Smoothed)' 

SAPQ13 = 'NEW CANCER MGT DELAY>14D S/P (Smoothed)' 

SAPQ14 = 'DIABETES UNCONTROLLED (Smoothed)' 

SAPQ15 = 'ADULT ASTHMA (Smoothed)' 

SAPQ16 = 'LOWER EXTREMITY AMPUTATION (Smoothed)' 



drop i; 



run; 



■*/ 
*/ 
■*/ 



/* 
/* 



merge rates 



data out3 . &outf ila3 

merge in3 . &inf ilea2 . (where= (_type_ in (0 , 8) ) 

keep=%area _type_ tapq01-tapgl6 papq01-papql6 oapq01-oapql6) 

riskadj (keep=%area rapq01-rapql6) 

smoothed (keep=%area sapq01-sapql6) ; 
by %area; 
run ; 



/* 
/* 
/* 



-*/ 
*/ 



output rates 
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proc means data=:out3 . fcoutf ila3 . (where= (_type_ in (8))) n nmiss min max 
mean sum; 

title4 'summary of area-level rates (_type_=8) ' ; 
run; 
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CANCER POP. DATA 



100 


35 


12083 


101 


20 


11697 


102 


15 


5046 


103 


30 


59218 


104 


20 


9505 


105 


25 


113419 


106 


20 


8669 


107 


90 


12083 


.108 


10 


1383 


109 


10 


8669 


110 


5 


14562 


111 


10 


9505 


112 


45 


8669 


113 


15 


12083 


114 


25 


8239 


150 


20 


113419 


151 


45 


59218 


152 


8 


5452 


153 


45 


1503 


154 


120 


25010 


155 


15 


11337 


156 


5 


6536 


157 


5 


9505 


■158 


5 


5163 


159 


25 


8669 


^160 


45 


78687 


161 


30 


6536 


162 


5 


10505 


163 


15 


4106 


164 


40 


12083 
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